function [tree, counts, scores, traits, exclScores, exclTraits, exclGenes, discScores, discTraits, singleScores, singleTraits, singleGenes, mixedScores, mixedTraits, almdiscScores, almdiscTraits, sharedScores, sharedTraits, sharedGenes, lowmixScores, lowmixTraits, subdiscScores, subdiscTraits, traitMatrix] = runAll(copyfilename, cancerfilename, subpops, maxGap)
%UNTITLED Summary of this function goes here
%   Detailed explanation goes here

    fprintf('reading in data...\n');
    [head, abs, rel, chr, chrEndSites, copy] = readData(copyfilename);
    counts = getCounts(subpops, head);
    fprintf('building transition vectors...\n');
    fprintf('finding events...\n');    
    chrEndSites = [0;chrEndSites];
    traits = [];
    scores = [];
    traitMatrix = [];
    for i = 2:size(chrEndSites, 1)
        fprintf('collecting traits from %d to %d\n', chrEndSites(i-1, 1) + 1, chrEndSites(i, 1));
        trans = getTransitionVectors(copy(:, chrEndSites(i-1, 1)+1:chrEndSites(i, 1)));
        [codes, starts, ends] = collectAllCombinedTraits(trans);
        fprintf('found %d %d %d traits\n', size(codes, 1), size(starts, 1), size(ends, 1));
        [traitsNew] = blendCombinedTraits(codes, starts, ends, maxGap);
        fprintf('found %dx%d unique traits\n', size(traitsNew, 1), size(traitsNew, 2));
        traitsNew = [traitsNew, ones(size(traitsNew, 1), 1)];
        for j = 1:size(traitsNew, 1)
            traitsNew(j, 6) = i - 1;
        end
        scoresNew = getSubpopulationCombinedTraits(subpops, trans, head, traitsNew);
        traitMatrixNew = getTraitMatrix(trans, traitsNew);
        
        traitMatrix = [traitMatrix, traitMatrixNew];
        scores = [scores; scoresNew];
        traits = [traits; traitsNew];
    end
    fprintf('found %d unique events.\n', size(traits, 1));
%     fprintf('calculating subpopulation scores...\n');
%     for i = 2:size(chrEndSites, 1)
%         trans = getTransitionVectors(copy(:, chrEndSites(i-1, 1)+1:chrEndSites(i, 1)));
%         thisChrTraitsIdx = 0;
%         for j = 1:size(traits, 1)
%             if (traits(j, 1) == i)
%                 thisChrTraitsIdx = thisChrTraitsIdx + 1;
%                 thisChrTraits(thisChrTraitsIdx, :) = traits(j, 2:size(traits, 2));
%             end
%         end
%         scoresNew = getSubpopulationCombinedTraits(subpops, trans, head, thisChrTraits);
%         scores = [scores; scoresNew];
%         
%         traitMatrixNew = getTraitMatrix(trans, thisChrTraits);
%         traitMatrix = [traitMatrix, traitMatrixNew];
%     end
    
    fprintf('Categorizing traits...\n');
    
    [exclScores, exclTraits] = getTraitsExclusive(scores, traits, counts);
    [discScores, discTraits] = getTraitsDiscerning(scores, traits, counts);
    [singleScores, singleTraits] = getTraitsSingleton(scores, traits, counts);
    
    [mixedScores, mixedTraits] = getTraitsMixed(scores, traits, counts);
    [almdiscScores, almdiscTraits] = getTraitsAlmostDiscerning(scores, traits, counts);
    [sharedScores, sharedTraits] = getTraitsShared(scores, traits, counts);
    [lowmixScores, lowmixTraits] = getTraitsLowMixed(scores, traits, counts);
    [subdiscScores, subdiscTraits] = getTraitsSubdiscerning(scores, traits, counts);
    
    [geneNames, geneChroms, geneStarts, geneEnds] = readCancerGenes(cancerfilename);
    
    if (size(exclScores, 1) < 1)
        exclGenes = [];
    end
    if (size(singleScores, 1) < 1)
        singleGenes = [];
    end
    if (size(sharedScores, 1) < 1)
        sharedGenes = [];
    end
    for i = 1:size(exclScores, 1)
        exclGenes{i, 1} = getGenesForTrait(geneNames, geneChroms, geneStarts, geneEnds, rel, exclTraits(i, :), chrEndSites);
    end
    for i = 1:size(singleScores, 1)
        singleGenes{i, 1} = getGenesForTrait(geneNames, geneChroms, geneStarts, geneEnds, rel, singleTraits(i, :), chrEndSites);
    end
    for i = 1:size(sharedScores, 1)
        sharedGenes{i, 1} = getGenesForTrait(geneNames, geneChroms, geneStarts, geneEnds, rel, sharedTraits(i, :), chrEndSites);
    end
    
    for i = 1:size(exclTraits, 1)
        chrom = exclTraits(i, 6);
        if (exclTraits(i, 2) == -1)
            startBP = rel(chrEndSites(chrom, 1) + 1, 1);
            startWin = 0;
        else
            startBP = rel(chrEndSites(chrom, 1) + exclTraits(i, 2));
            winTop = rel(chrEndSites(chrom, 1) + exclTraits(i, 2) + exclTraits(i, 3));
            winBot = rel(chrEndSites(chrom, 1) + exclTraits(i, 2) - exclTraits(i, 3));
            startWin = max(winTop - startBP, startBP - winBot);
        end
        if (exclTraits(i, 4) == -1)
            stopBP = rel(chrEndSites(chrom + 1, 1), 1);
            stopWin = 0;
        else
            stopBP = rel(chrEndSites(chrom, 1) + exclTraits(i, 4));
            winTop = rel(chrEndSites(chrom, 1) + exclTraits(i, 4) + exclTraits(i, 5));
            winBot = rel(chrEndSites(chrom, 1) + exclTraits(i, 4) - exclTraits(i, 5));
            stopWin = max(winTop - stopBP, stopBP - winBot);
        end
        exclTraits(i, 2) = startBP;
        exclTraits(i, 3) = startWin;
        exclTraits(i, 4) = stopBP;
        exclTraits(i, 5) = stopWin;
    end
    
    for i = 1:size(sharedTraits, 1)
        chrom = sharedTraits(i, 6);
        if (sharedTraits(i, 2) == -1)
            startBP = rel(chrEndSites(chrom, 1) + 1, 1);
            startWin = 0;
        else
            startBP = rel(chrEndSites(chrom, 1) + sharedTraits(i, 2));
            winTop = rel(chrEndSites(chrom, 1) + sharedTraits(i, 2) + sharedTraits(i, 3));
            winBot = rel(chrEndSites(chrom, 1) + sharedTraits(i, 2) - sharedTraits(i, 3));
            startWin = max(winTop - startBP, startBP - winBot);
        end
        if (sharedTraits(i, 4) == -1)
            stopBP = rel(chrEndSites(chrom + 1, 1), 1);
            stopWin = 0;
        else
            stopBP = rel(chrEndSites(chrom, 1) + sharedTraits(i, 4));
            winTop = rel(chrEndSites(chrom, 1) + sharedTraits(i, 4) + sharedTraits(i, 5));
            winBot = rel(chrEndSites(chrom, 1) + sharedTraits(i, 4) - sharedTraits(i, 5));
            stopWin = max(winTop - stopBP, stopBP - winBot);
        end
        sharedTraits(i, 2) = startBP;
        sharedTraits(i, 3) = startWin;
        sharedTraits(i, 4) = stopBP;
        sharedTraits(i, 5) = stopWin;
    end

    fprintf('Generating Tree...\n');
 
    tree = generateTree(traitMatrix, head);
end

function genes = getGenesForTrait(geneNames, geneChroms, geneStarts, geneEnds, rel, trait, chrEndSites)

    chrom = trait(1, 6);
    if (trait(1, 2) == -1)
        startBP = rel(chrEndSites(chrom, 1) + 1, 1);
    else
        startBP = rel(chrEndSites(chrom, 1) + trait(1, 2) - trait(1, 3));
    end
    if (trait(1, 4) == -1)
        stopBP = rel(chrEndSites(chrom + 1, 1), 1);
    else
        stopBP = rel(chrEndSites(chrom, 1) + trait(1, 4) + trait(1, 5));
    end
    
    genes = '';
    for j = 1:size(geneNames, 1)
        if (geneChroms(j, 1) == chrom && geneStarts(j, 1) > startBP && geneEnds(j, 1) < stopBP)
            genes = strcat(genes, geneNames{j, 1}, ', ');
        end
    end
end

